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Abstract 

The methods of statistical dynamics are apphed to a fluid with 5 
conserved fields (the mass, the energy, and three components of mo- 
mentum) moving in a given external field. When the field is zero, 
we recover a previously derived system of parabolic partial differential 
equations, called 'corrections to fluid dynamics'. 

Key Words: Navier-Stokes, bulk diffusion, Soret, Dufour, thermal diffu- 
sion. 



1 Introduction 



In a previous paper |23|, we described a stochastic model of a fluid with 
no external field and derived a system of parabolic equations expressing 
the dynamics of the density-fields of mass, energy and momentum. An 
unusual feature of the system is that the Euler continuity equation acquires 
a diffusion term. This bulk diffusion does not appear in the standard theory 
, but has arisen in some other work ^, ^ . It is likely that our equations 
are more stable and physically more accurate than the usual Navier-Stokes 
equations. In a model in an external field, but without a velocity field |21], 
we studied the dynamics of the mass and energy densities; they differed 
from the equations got by putting the velocity field, u(x, t) equal to zero 
in the Navier-Stokes system. In particular, the mass-density p{x, t) did not 
obey the Euler continuity equation with u = 0, but rather the Smoluchowski 
equation. This contains not only the diffusion term found in [23|, but also 



a drift, as predicted by Smoluchowski, and implicit in the work of Einstein. 
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In the present paper, we take the model of 
field We derive the full set of equations 



|, and put it in an external 



I + div(u,) = A div (p-iv(evv) + ^ 



(1) 



djpe) 
dt 



dpu 
~dt 



+ div [\i{pe + P)] = \ div \2p~^V [q^I'^p] + p~^V [Q^'^p] 



+ 



+ 2—-r^P 



+ div(pu (g) u) = -pV<I>/m - VP + A div 



u 



/2, 



(2) 



(3) 



Here, p is the mass density, u is the velocity field, and is the temperature. 
The energy density per unit mass, e, reduces in this approximation to <I>/m+ 
3fcs0/(2m). The equations reduce to those of |23| when <I> = 0, while 
if u = they reduce to with suitable changes due to the different 
multiplicity of states in the two models. This is a property of any good 
approximation method, which is not true of the usual compressible Navier- 
Stokes equations with temperature, as studied in [l^, 17 1 for example. The 
equations of [^] can be extended to models with interparticle forces, treated 
as a nonlinear mean field theory by a development of the theory of non-linear 
parabolic systems |Q, ||, |6|, |^, |^ . 



2 Review of the model 

2.1 The sample space, conserved quantities and information 
manifold 

Let A be a finite subset of the cubic lattice (aZ)^ with spacing a ~ 10~^ cm, 
representing the size of the hard core of the fluid molecules. We endow our 
model with a local structure by assigning to each x £ A a sample space 
which specifies the possible configurations at x. The total sample space is 
then the product space il = YlxeA^x- For a structureless monatomic fluid 
(e.g. argon) we choose 

Jl, = {0,(eZ)3}, (4) 

where e ~ 6.6 x 10^^^ c.g.s. is the quantum of momentum of a particle 
conflned to a region of size (see the discussion following equation (6) in 
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[p3[) . A point u; G 17 is thus the collection {uJx}x£A- If ^ is such that uJx = 9 
for a certain x £ A then the configuration u! has no particle at x; we also 
say that there is a hole, or vacancy, at x. If w^,' = k S (e^)^, then in the 
configuration to there is a particle at x and its momentum is k. 

The state space of the system consists of the probability measures on Q 
and is denoted by The information manifold 5 C S associated with 

the model consists of the parametric exponential family determined 
by the slow variables, that is, the local mean values of the five conserved 
quantities: the mass, energy and the three components of the momentum 
along the three unit vectors 61,62,63 generating the lattice. These are the 
random variables 

^-^"^ = I 1 if... = kE(6Z)3 

r. / X _ / iiujx = 9 

■^^^""^ ~ \ (2m)~ik.k + $(2;) ifu;^ = k 

^ I . / if = 
= jk if.. = k, 

where m is the mass of the molecule and <I> is a given real-valued potential, 
which could be time-dependent. The states in 5 C S are those of the form 
P = WxPx^ where 



Px 



exp { -ixMx - Px£x -Cx-Tx]. (8) 



Here, S,x,Px and Cx are fields of intensive variables that obviously determine 
p; they are called the canonical coordinates for p £ S. The great grand 
partition function Sx at each site is the normalising factor 

Ex = Yl ^^P{-^xMx- PxSx-Cx-'Px} (9) 



1 + e-«--^-*(^)ZiZ2Z3, (10) 



where 



This can be calculated explicitly if, due to the small parameter e, we ap- 
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proximate sums by integrals: 

e-%r-afeci/c (12) 

-oo 



TJ 



1/2 



^m{cr/m_ (13) 



The means {Nx,Ex,vJx) of the slow variables in a state p e M also 
determine p; these are called the 'mixture' coordinates. They are related to 
the canonical coordinates by a Legendre transform: 

= xGA (14) 

K = ^ = 1,2,3, xeA. (16) 



Using the explicit expression for the partition function, these can be used 
to deduce several equations relating the macroscopic variables in the theory. 
In particular, if we introduce the mean velocity field 

and the temperature field 6a; = l/{kBl3x), it is straightforward to show that 
and 

= (^^x) + ^kB@x + ^mux ■ vLx^ . (19) 

The discrete nature of the model allows us to use the von Neumann 
entropy 

S{p) := -kB ^p(u;) logp(a;), (20) 

which for p ^ S gives 

S = Y1 (^xA^x + /?x^x + Ca= • -C^a, + log Sx) (21) 

xeA 
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An argument from equilibrium theory then leads to the definition of the 
thermodynamical pressure as 

P(x) = 0-3^:56 log Hx. (22) 

If the ratio Vq/V (between the smallest volume that the = Nx parti- 
cles can occupy, that is Vq = a^N, and the total volume V) is small, then 
the formula above for the pressure reduces to the perfect gas approximation 

P. = (23) 

The other macroscopic variables in terms of which the hydrodynamical 
equations are written are the mass density p{x) = mNx/a? and the energy- 
density per unit of mass e{x) = Ex/{mNx). If we ignore the small term 
involving u^; • u^; in (|l^), we see that 

2.2 The hopping dynamics and the continuum hmit 

Whether time is discrete or continuous, the dynamics traces out an orbit 
in S. The need to take the time- interval larger than zero has been well 
explained in For discrete time, the dynamics consists of two steps. The 
first one is a stochastic map T : S — > S which maps each shell 

^N,E,z:7 =<uj€n: = N,^£x=E,^Vx = ^\ 

\ X X X ) 

to itself. This is obtained when we specify the hopping rules that are re- 
sponsible for changing the configuration of the particles and hence the local 
values of the slow variables. If we start with a state p £ S, the state T{p) 
will generally not be an exponential state. The second step in the dynamics 
is then a orthogonal projection back to S following a path that conserves 
the means of all the slow variables. This is the thermalising map Q, de- 
fined for any state p with finite expectations of the five slow fields, as the 
unique point in S with these expectations as mixture coordinates. So Qp is 
determined by 

(25) 
(26) 
(27) 



^QpWx] 


= ^pWx] 


= Nx 


'^Qpi^x] 


= Ep[<5x.] 


= Ex 


^Qpi'Px] 


= ^p[Vx] 


= '^x 
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The precise meaning of orthogonality in the description above, as well as the 
characterisation of the path followed during the projection as a geodesic for 
a certain affine connection (the 'mixture' connection), is part of the subject 
called information geometry |§, 12, |l3|, |l5|. 



Our model is specified by giving hopping rules. We require that T should 
couple only neighbouring points in A, where we consider two points to be 
neighbours if their distance along one of the lattice unit vectors is one mean 
free path, denoted by i. A particle moving from a site is taken to move 
exactly the distance i, and then to thermalise. This was called abrupt ther- 
malisation. The more elaborate assumption, that the size of its hop is ran- 



dom, governed by the exponential law with mean £, leads [24| to similar 
conclusions. We assume that i is an integer multiple of the lattice spacing 
a, but allow it to depend on the local density by taking i/ a to he the nearest 
integer to 



m 



p{x) a^p{x) 

Suppose that a; G 17 is such that x G A is occupied. Consider in turn the 
possibility of jumping from x along the direction of the unit vectors itej of 
the cubic lattice to an empty site x' := xzizlci, i = 1,2, 3. In the absence of an 
external potential |2^], the jump will take a time where Vj, := k^^/m. 

We then define the (random) hopping rate from x to x + iei to be the inverse 
of this relaxation time, namely vl./i if wj. > and zero if is negative, in 
which case there is a rate —v\li of hopping to x — £ej. The situation in 
the presence of an external potential is a little more involved, because the 
potential causes a change in the velocity along the jump. From now on, 
assume for definiteness that <^(x + le,-i) > ^{x). If the particle at x hops to 
x + iei, its potential energy increases to ^{x + iei) and so its kinetic energy 
must decrease by the same amount. Its change in momentum is taken to be 
entirely in the direction of e^. So its velocity in the i-direction, Vx ■= kx/m 
(we omit the index i for v and k in the following formulae since it is clear 
from the rest of the notation what is the component involved) is reduced to 
v'x := k'x/m, where 

k'^/{2m) = kll{2m) - ^{x + iei) + $(x) 

= kl/{2m)-ldMx) + 0{e), (28) 



that is 



k'^ = {kl - 2midMx)Y''^ . (29) 
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In order for the move to be energetically possible, we must have 

fcx > 4 := {2imd^^{x)f/^ . (30) 

Similarly, if the particle at x with velocity Vx < in the i-direction hops to 
X — ici, its potential energy decreases to <^(x — isi), with a corresponding 
rise in its kinetic energy. Therefore (again taking the change in momentum 
to be entirely in the i-direction) its (negative) velocity in the i-direction 
becomes v" = k"/m, where 

A:f /(2m) = kl/{2m) + <^{x) - $(x - £e^) 

= kl/{2m)+£di^{x-eei) + 0{f), (31) 



that is 



k'^ = -[kl + 2m£di^{x - ici)) ^''^ . (32) 



We take the hopping rate from x to x' = x zt (.Ci to be the average of the 
initial and final rates: 



r{kx 



r+[Kj:) .— 2mi ~ 2mi ' ^ Kx ^ K^. 



These hopping rates increase with kx and there are infinitely many possible 
momentum states. To be a Markov chain, the sum of all rates out of a 
configuration must be less than one. For any kx this can be achieved by 
choosing dt small enough. To do this for all kx with a fixed dt we must 
put in a cut-off; there are no hops if \kx\ > Kx, say. Finally, r{kx)dt gives 
the probability of a transition in an interval dt provided that the site x is 
occupied and the site x' is empty, so the actual entries of the Markov matrix 
are conditional probabilities and the transition rate above should appear 
multiplied by factors of the form A'^j,(l — Nx'). As argued in |^, we neglect 
Nx comparing to 1, therefore leaving out the second term in the factors 
above. 

The continuum limit we are going to take in order to obtain the hydro- 
dynamical equations corresponds to £ — > 0, c — oo such that the product £c 
remains finite and non-zero, where 

.:=fM£y'' (34) 

\ m J 
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is the approximate velocity of sound at the reference temperature Qq. The 
diffusion constant that appears when we take the hmit is then predicted to 
be 

^ icp ^ aPmaxC , . 

(27rGo)i/2 (27rGo)i/2- ^""^^ 

3 Hydrodynamics in an external field 

When a transition from x to x + ici occurs in a potential <1>, the loss of 
mass and energy from the site x is equal to the gain at the site x + ici. 
This is not true of momentum; the loss at x differs from the gain at x + ici 
by Ki := {2£mdi^)^^'^ . So we deal will M and £ first. We take it that if 
< k' < Ki then no hop is made. 

Before we start the calculations, let us recall that integrals of the form 

POO 

M„(C) = / k''exp{-pk'^/{2m)-Ck}dk, n = 0,l,2,3 
Jo 



were evaluated up to second order in C, in Appendix 1 of [g3| . For later use, 
we reproduce the results here up to zeroth order in for n = 0,1, 

M„,C) ^ {^f, (36, 

777 

MiiC) = J, (37) 
and to first order in for n = 2, 

3.1 Dynamics of the mass-density in an external field 

Since the field <I>(x) is external, it does not depend on the configuration oj^ 
of the random fields at x. It therefore cancels in the exponential states. 
This is seen in its simplest case in the model studied in [^]. The potential 
enters only in its supression or enhancement of the transition rate; in the 
present case, the rate is the average of the initial and final rates. This shows 
up mainly in the appearance of a non-zero lower limit to the (positive) 
momentum for any right-going hop to be possible. 



8 



Let Jx/l be the change in the value of due to the hoppings occuring 
between x and x — ici in such a way that the change due to exchanges with 
both X lb £ei in an interval 6t is 

S^N^ = -^^^±^^p^5t. (39) 

So the total change in A'^,' in an interval 5t due to hoppings in all directions 
is 

6N, = (diN, + 62N, + 63N,)6t. (40) 

Using the hopping rates defined in the previous section, the loss/gain 
contribution to the particle current involving the exchange between x and 

X — iei is 

x—lei 

(41) 



As in |2^, the analysis of this expression is best handled by introducing 
a conditional probability PxCf-^) = V^^\^x = 1) on the particle space 17 — 0, 
that is, 

p,(k) = (ZiZaZg)-! exp{-/3,|k|V(2m) - ^ • k}. (42) 

We now use the fact that Px{k) = Nxpx{k) and Nx{k) = 1 on the particle 
space — 0, replace the sums by integrals in (pl), and add and subtract the 
term 

Fi = NxiZie)-' j lr+{k) exp - C^j dk, (43) 

to obtain 

We start by calculating F^. In the term 

k>K 2m [ 2m J 
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we make the change of variable 



1 12 1 2 2 

k = k — K . 



Then kdk = k' dk', and the integral becomes 



exp 



,>o 2m(A;'2 + K2)i/2 

which can be written as 

k' [(k'^ + ^2)1/2 ^ j^r 



^(^_±^_^^(fc'2 + ^2)l/2|^^/^ 



k>>o 2m(fc'2 + k2)1/2 



exp 



pk 



12 



2m 



Ck' 



X exp 



2m 



(^(fc'2 ^ ^2)1/2 _ ^/^ ^^ I 



The arguments of the exponentials are small, and we expand them to first 
order: 



exp 



2m 



fA;'2 + K2^1/2_^, 



1 



This gives us the three terms 

k' [(A;'2 + ^2)1/2 ^ ^/j 



exp 



fc,>o 2m{k'^ + ^2)1/2 

/■ fc' + At2)l/2 + ^/ 



2m yfc,>o 2m(A;'2 + k2)1/2 

exp 



2m 



2m 

exp 
/9;^'2 



A;' - ik'^ + 



/2 , ,.2^1/2 



(45) 



13k 



12 



2m 



2m 
Qk'X dk'. 



Ck' } dk' (46) 
(47) 



2m Jk'>o {k'^ + ^2)1/2 
The dominant term is (|45|), in which we may replace the factor 

k' [{k'^ + ^2)1/2 + 
2m(fc'2 + ^2)1/2 

by the velocity when $ = 0, namely k'/m, with an error of 0(^log£) (Ap- 
pendix 1). So the contribution of this term to the mass current can be 
approximated in the limit by 



Zie Jo m 



(3k 



12 



2m 



Ck' \ dk. 



(48) 
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Making the same replacement in (|4q), with the same error, we obtain 

' — exp<^-^- Ckfdk 



2m Jq m I 2m 



2m^ 

Therefore, the contribution coming from ( ^6|) to the mass current is 



which, to zeroth order in gives 

1 /2 

_Ar^ ( JL\ idMx) = ^-^dMx) = 

(49) 

When multiphed by mjc? this is what we call the Smoluchowski, or drift, 
current: 



Mo(C). 



The integral in ( [47|) is bounded by 



2m 

Expanding Mq to zeroth order in Q ^ this gives 



2m 2^/37 2 V/3y V2y c 

so that it can be ignored in the limit. 

We now turn our attention to the second term in the current ( ^4[ ) . We can 
again replace the factor by |A;*|/m (Appendix 2), so that the contribution 
to the mass current from this term is 

f r ^e.p{-f!-ca*. (51) 

■^iC y_oo m y 2m J 

When we combine this with (45) and multiply it all by mjc? what we find 
is simply piU' . 



— 9,«>, 
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Finally, we need to deal with the last term in (44), which in the limit 
becomes just —idiF^. As we have just shown, the only non-negligible terms 
in itself are ( ^8[ ) and (|4^) . But ( ^9| ) is already of order £c and therefore 
can be ignored when multiplied by the additional £ above. The only term 
that survives is 



-idi 



Zie 



k 



m 



■ exp 



2m 



Ck } dk 



m V 27rm 
which, to zeroth order in (^*, is 



1/2 



di 



( (2^/3)1/2 ; (27reo)i/2 y ^1/2^1/2 



-^di{Q^'^N^). (52) 



When multiplied by mjc? this is what we call the diffusion current 



A 



5i(ei/V), 



which is made up of the Fick current 

-Ae^/2^,(iogp) 

and the Soret current 

-A(2G^/2)"^9ie. 
Therefore, we obtain the total mass current by collecting together 



^ %9.cl>(x) - ^5.(ei/2iV.). 



(53) 

(54) 
(55) 



(|lD, (gD and (|5|), that is 



pkeS^/^ 

We now go back to (p9|) and expand the finite difference in there as 



(56) 



p _ p 



dJ, 



+ 0{t). 



Since the expression we just found for J* does not contain any term with a 
large factor c, we see that, in the limit £ — > subject to keeping ic finite, 
equation (^) becomes 



dt 



+ divJ = 0. 



(57) 
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Multiplying both sides of (|^) by m/a'^ gives us the equation for the time 
evolution of the particle's density 

1^ + div(up) = A div (p-iV(eVV) + , (58) 

or 

^ + div( Jp) = (59) 
where the conserved density current is found to be 

Jp = up+ Jd + Js- (60) 

3.2 Dynamics of the energy in an external potential 

Let Jx/^ be now the change in the value of Ex due to the hoppings occuring 
between x and x — iei. As before, the change due to exchanges with both 
X lb ici in an interval 5t is 

S^Ex = -:^^±^|— ^<5t (61) 

and the total change due in in an interval 6t due to hoppings in all 
directions is 

SEx = {SiEx + 62EX + 5^Ex)5t. (62) 

We have that 

4 = - ^ tr.{kx)p:,{y)Sx{y) + er+{kx-eeJp.-eeA^)£.-eeA^), (63) 

where f = k • k/2m + ^(x). 

The analogue of the quantity Fx of the previous section is now 

Gi = Nx{Ze^)-^ I £r+{k') + ^{x)] exp S-^^^^ - C • kl d^k. 



2m I I 2m 



(64) 



Adding and subtracting this to (p3|), replacing sums by integrals and 
again using that Px{k) = NxPxik), we obtain 

-i{Gl-Gi_,,J/t (65) 
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We calculate Gx first (for i = \). In the integral 



ki>K 



/ti + ft? - k2)V2 /k-k , A 



2m 



/3k k 
2m 



C • k ^ d^k, 



we make the change of variables A;^^ = k\ — k? while keeping ^2 = ^2 and 



k'^ = k^. Note that 



and kidki = k'ldk'i, so defining 



A(k') 



k^-k^ 

2m 



+ ^(x + £ei) ) exp 



/3k' • k' 
2m 



the integral becomes 
,;>0 2m(fcf + «:2)i/2 



A(k') exp 



2m 



(66) 



We now expand the exponential to first order 



exp 



2m 



2m 



Thus the integral splits in the following three terms 



r k[ [(fcf + K2)V2 + fc;] 

Ar,/3K2 /■ k[[{k[^ + K^)^/^ + k[ 



Ze-> A,>o2m(fcf + k2)1/2 



A{kl)d^kl. 



(67) 
(68) 
(69) 



If, we approximate the hopping rates appearing above simply by k'^/m 
(Appendix 3) and use that <I>(x + ^ei) = $(x) +^9i$(x), the integral in ( |67| ) 
becomes 



Nx 


1 ^ 






Nx 


[ K 


Ze3 


Jk[>o rn 



k[ fk' -k' 



m V 2m 



+ ^{x)j exp 

Pk' ■ k' 

2m 



(3k' ■ k' 



2m 
C-k'\ d^k' 



C-k'\ d^k' 



(70) 
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The first part of (70) is later going to be combined with the integral over 
the negative values of ki appearing in the second term of (|65|). As for the 
second part of (|70|), we have 



f k'l 0^ ^ 



k[>0 



/3k' • k' 



Ziem 



2m 

12 



C • k' [> d^k' 



Ziem m 
Expanding it to zeroth order in C^, we get 



27rm 



Mi(C^)e ^ 



NJdMx) 



27rm 7 



-N^Q^/^di^x). 
P 



(71) 



For the integral (^), we again approximate the hopping rate by k'^/m 
and use that $(x + ^ei) = <I)(x) + idi<^{x), to obtain 



k[ 



k' k' 

2m 



exp 



/3k' • k' 
2m 



C • k' !> d^k' 



/" A;'i($(x) + ldi^{x)) exp | -^^-^ - C • k'l (i3k^72) 
Jk[>o I 2m J 



The first term above divides into three integrals. The first one is 



2m'^Ze^ 



2m 



exp 



/3k' • k' 
2m 



C • k' !> d'k! 



the second one is 



mZie 2m 
/3 V^' 



27rm 



2m2ZiZ2e2 



£di<^{x)kB@; 
kf _^Aifc, „ 



(73) 



k[>0 



2m 



0k 



kle-^-^^-'dki 



k2 



mZ\e 2m 

1/2 



m^{C 



2'Km) ^ ' 



/32 

m{v?)'^ 



m 

^1. 



1 



+ -ksQ 
2 2 



(74) 
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while the third one is 



2m?ZiZ-ie'^ 



-J—e 2m '^1 



k[>0 



2m 



dk'\ 



mZi€ 2m 



m2(C3)2 



/32 



m 

+ 7 



1 



2 +2^-® 



2'irm 

So the total contribution from the first term of { \l2\ j is 



27rm 



(75) 



(76) 



The terms involving the velocities disappear in the limit. The remaining 
term is 



-2-N^Q^I'^di^{x). 
P 



As for the second part of ( [72D we have 



2m'^Ze^ 



Mx)+£di<S>{x)]Mi{C^) 



2TTm 



1/2 



(77) 

idi^{x) [$(x) + £9i$(x)] 
(78) 



of which only the first term survives in the limit, leaving us with 

The integral in ( |69| ) is itself of smaller order and can be ignored (Ap- 
pendix 4). This completes the contribution of G]. to the energy current. 

We move on to deal with the second term in (|65|). The factor £r_ can 
be replaced by —ki/m (Appendix 5), leading to a contribution of the form 



Ze3 



k 



fci<0 



kk 



m \ 2m 



+ <l>(x) I exp 



^k-k 

2m 



As promised, this joins the first part of (|70| ) to give 



m 



m 



C • k ^ d^k. 



pi 



(80) 



m 



(81) 
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We finally look at the last part of (65), which in the limit becomes 
—idiGx- As we have just seen, all the contribution coming from Gx are 
already of order ic (and can therefore be discarded when multiplied by the 
additional £ above) with the exception of the first part of (|70|), that is, 



Nx_ 
Ze3 



k[ /k'-k' 



k[>0 



2m 



+ ^{x) exp 



/3k' • k' 
2m 



(82) 



We recognise the term not involving $ as being equal to the first term in 
( [72D times a factor —2m/ (3k? = —[j3i.di(^{x))~^. Therefore, its contribution 
to the energy current is 



-Idi 



A 



2^NxkBe>'i' 



-2-di{NxkBQ'"^). 
P 



(83) 



As for the term involving <I>, we again recognise it as being equal to the 
term involving <^ in the second part of (^) times the factor —2ml (3k? = 
— l//3£0i$. Therefore, its contribution to the energy current is 



--5i(iV^ei/2)$(x) _ -f^^Q^I'^dx^{x) (84) 

p p 



We can now take a breath and collect all the terms we have obtained for 
the energy current to put back in (||). They are (0), (0), (|7|), (|8l|), (g) 
and (|8^ ) and the end result is 

Jl = u'x{Ex + NxkBQx)-2-^NxQ^'^d^^x) 

- ^TtM— ^(^) - '^-m.kBQ'") - -d.{Nxe'/'mx). (85) 
kB^ ' P P P 

Again we see that none of the terms in J* contains a large factor c, so 
that the discussion preceding (^) applies here as well and in the limit we 
obtain 



dt 



+ divJ = 0, 



(86) 



that is 



dEx 

dt 

V^(x) N. 



+ div [ux{Ex + NxkBQx)] = div 



2-Nxe^^^V^x) 
. P 

+ ^T^T^— ^(^) + 2-V(iV,A;Be3/2) + -V{Nxe'/^)Hx) 
kBK^ ' P P P 



(87) 
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We now use that e{x) = Ex/{mNx), 4>{x) = $(x)/m, p{x) = mNx/a? 
and Px = NxKb@x/cl'^ , divide both sides of the previous equation by and 
obtain 
d{pe) 



dt 



+ div [u{pe + P)] = X div 2p-^V (^B^/^p^ ^ ^-ly 



which can be written as (recall ( |50D and (|53|)) 

f pel/2 



+ div [u(pe + P) + {Jd + Js)0 + 2J5P] = 2A div 



(89) 



3.3 Dynamics of the momentum in an external field 



Since momentum is not conserved (as there are body-forces due to the exter- 
nal field), the rate of change of momentum density will not be the divergence 
of something; we expect the extra term to be pf where f := —V^/m is the 
force per unit mass. To see this, let us define the current J* as in the 
previous sections, namely 

j; = - £r_{kMk)Vi{k) + ir+{kl^,Jpx-ee,{k)vi_,Jk) (90) 

Then the change in wi due to exchanges with both x it ici in an interval 
5t will only be given by the usual 

S^^J^ = J.kM_Ast (91) 

for i ^ because it is implicit in this formula that the particles hopping 
from X lb ^e, to x have their j-component of the momentum unchanged 
during the jump. We do this case first. The analogue of PJ and G]. from 
the previous sections is now 

Hi = Nx{Ze^)-^ [ er+{k')k^expl-^^^-C-'k\d^k. (92) 

Adding and subtracting this to (^) , replacing sums by integral and using 
that Px{k) = NxPxik), we obtain the familiar form 

4 = Hl-^l ^r_(F)fe^exp|-^^-C-kl(i3k 



HI - Hi_,^J IL (93) 
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Notice that since Wx = mN^Ux, in this section we shall keep all terms of 
order mux- If we now perform in (p^) the integrations over kr {r i, r j) 
and kj we find 



jNx 
mui. 



k>K 



Ak) exp I —— C,k\ dk 

2m 



But now we can use the calculation we have already done for F*, that is, in 
the limit we have 



Hi. 



mu: 



Nx_ 
XNx 



°° k' 
— exp 
m 



f3k 



12 



2m 



Ck' } dk 



pkB@^/^ 



mui 



(94) 



The second term above, when multiplied by m/a^, simply gives muxJg. 

As for the second term in (p3D, we can also perform the integrations over 
kr (r 7^ i,r ^ j) and kj, as well as to replace the factor ir^{k'^) by —k^/m 
to find 



Nx 



k 



m 



f3k'^ 



2m 



(95) 



When we add this to what we have found in the first term in (^) and 
multiply them m/a?, the resulting term is mpuxu\.. 

Finally, we look at the last term in (^), which in the limit becomes 
—idiHl. Since the second term in (p^) is already of order £c, we see that 
the only surviving contribution here is 



mui 



k' 



m 



■ exp 



(3k 



12 



2m 



Ck' } dk 



-id. 



mui 



Nx 
' Ziem 



(96) 



which, to zeroth order in C gives 



Nx 



mui 



1/2 



m \ 2-nm , 



di{NxQ^'^mu^) 
d,{NxQ^''^mu^). (97) 



(2^eo)i/2 
A 
P 
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Therefore, the momentum current, for i ^ j, is given by 



J* = mNxul.ul. 



mui - -di{N^e^/^mu^). (98) 
P 



Once more, none of the above terms contain a large factor c, so in the hmit 
£ — > subject to ic finite, we can approximate the term (J*_,_^g. — Jx)/(^ in 
(|l|) simply by diJ^. 

For i = j, equation (|9T|) is not correct. In this case, if kx-ee^ > i^x-ie- is 
the j'-component of the momentum for a particle at x — iej, then 

k'x_ee^ = [kl_,,^ - 2mld,^{x - ^e,)] 

is the j-component of its momentum when it arrives at x. Similarly, if 
kx+ie < is the j-component of the momentum for a particle at x + 



then 

k'Uie,=-[kUie, + 2mid^^{x)]y^ 

is the j-component of its momentum when it gets to x. Therefore, the total 
rate of change in due to exchanges with x it Icj consists of the following 
four terms 



5jwi 
6t 



r-ikx+ie,)kx+ee,Px+ee,{K+ee,) (99) 
- ^ r+{kx)kxPx{^x) - ^ r_{kx)kxPx{^x) (100) 

+ Yl r+{kx-lej)k'^^f,^^Px~ie,iK-le,) (101) 



If we now recall from ( ps] ) what the hopping rates look like, and use that 

kx+ie, + ^"+^e^) = (kx+ie, + k'x+te^ kx+ie, + 2midj^{x) 

^ " + k'x_ie?\ K-ie, = [kx-ie, + A^^-^eJ kx-ie, - 2lmdj^{x - ICj), 
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then we can rewrite the above as 

"^Jf' = XI f-{^x+eej)kx+eejPx+eej{^x+eej) 

- X r+{kx)kxPx{K) - X r-{kx)kxPx{K) 

+ X] r+{kx-eej)kx-eejPx-eej{K-eej) 

- dj^{x) ^ Px+eej{kx+eej) 

kx+eej <0 



which we then recognise as 

5jwi _ Ji+eei ~ 



St e 



dj^{x-lej) Px-iej{K-ie,)- (102) 



fe ■ >K„ 



In the limit, the last two terms above add up to Nxdj^{x) (Appendix 6). 
Let us now calculate J|, that is 

4 = Hi-^J £r_(fc^)Fexp|-^-C-klci3k 



^ Ze^ JkJ<o { 2m 



^ {Hi - Hi_,^) 11. (103) 



We find that Ri reduces to 



^ Z,ey,>^ 2m ^1 2m ^ j ' 
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where we can make the famihar change of variables /c'^ = k'^ — k'^ to obtain 



Jk'>0 



(A;'2 + ^2)1/2 + ^/ 



2m 



k' exp 



pk 



n 



2m 



C^k' 



X exp 



2m 



{[k'^ + n'f/^ -k')Adk' 



Expanding the exponential to first order gives the usual three terms 



Hi 



2m 
' 2m 



k'>0 



k'>0 



{k'2 + ^2)1/2 ^ fe/- 

2m 

■(A:'2 + ^2)1/2 + ^/ 



Jk'>0 



exp 



2m 

_(3k^ 
2m 



k' exp 
k' exp 



Pk 



/2 



2m 

(5k'^ 
2m 



C^k' } dk' 



C^A:'|dA(l04) 

-C^A;'|(i/('l05) 
(106) 



In ( PI , we can replace the factor 

(A:'2 + ^2)1/2^^/ 



2m 



by fc'/m (Appendix 7), so that it amounts to 



(/fcO'exp 



f3k 



12 



In ( [L05[ ), the same replacement gives 



2m 



C^k' } dk'. 



2Zjem'^ 



k'>0 



{k'Yexp 



2m 



Ck' } dk' 



2Zjerpfi 



(107) 



M2(C^'), (108) 



which, to first order in (recall that we are keeping terms proportional to 
mv? in this section, since tz7^ = N^mux)., gives 



2Nx(3mldji^{x) ( P 



2m2 



2'Km 



2Nx 



P 



2TTm 



1/2 



7r'^i/2 ^ m\ 



j^[x)mu^ 



3/2 



21^1 



of which only the second term survives in the limit, resulting in 



-2 



Ql/2 



dj^{x) 



mu-i 



(109) 



(110) 
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Similarly in ( |106| ), we obtain 



2Zjem 



2m V 27rm 



1/2 



2iTm 



£dj^{x)mul. 



mui 



(111) 



Moving to the second term in ( |103| ), we replace the factor ir-{k^) by 
-h? jm (Appendix 8), so it contributes with 



Combining ( 107 ) and (|112| ), what we obtain is 



mZje 



2m m. 



i\2i 



m 

= mN^iuif + N^kee^. (113) 

Finally, for the last term in ( |103D , which in the limit becomes —idjHx, 
only (107) contributes, since the other terms in are already of order ic. 
We get 



-idjHi 



-id, 



Zjcm 



k'>0 



{k')"^ exp 



2m 



- C^k' } dk' 



Zjtm 



which, to first order in (^^ , gives 



1 m \ 2'iTm 



1/2 



7r\ 1/2 f m\ 



3/2 / N 2 ■ 



'-dj{N^kB@^)-2£dj 



mu 



P \ 1 
27rm ) (3 



(114) 
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Therefore, the momentum current, for i = j is given by the sum of ( pT 
(Iml ), ( pl3|) and (|ll|), that is 



dj^{x) 



mut 



We see that the term N^ksOx above is of a larger order than the others 
and is not neghgible in the expansion of the finite difference ( Jx+£e ■ — Jx) 1^- 
What we obtain is 



a jj p pp. jj 



Thus, in this case, the finite difference (^Jx+ee — Jx) can be approximated 
by 



djlmNxiui)^ + NxkBQa 



A 



mu: 



mut. 



(115) 



Therefore, collecting together the contributions from i = 1,2,3, we see 
from (^), ( |102[ ) and the equation above, that the change in vj^ is governed 
by the equation 



d 



dt dxi 

+ NxkB@x^ij 



pdx^ \ 



pkB@^/^ J pdx 



d 

Nxdj^{x) -Y^Q-i 
A d 



ut 



i=l 



(116) 



The first term above is not covariant and we need to average it over 
50(3). The averaging procedure is explained in |23|, and its result is 



V— — 

^ dx^ bp 

1=1 



dx^ 



— N.xQ'/'mui)+2—NxQ' 



d 



dx^ 



1/2 



(117) 



We can now put it back into ( |116| ), divide both sides of it by and 
use that vui = mNxui, p{x) = mNx/a?, fi = —dj^{x)/m and Px = 
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NxkB&x/O'^ to obtain 

dp{x)ui 



dt 



2 



In vector notation, this reads 
dpu 



(118) 



dt 



+ div(/9U0u) = pf- VP + Adiv( ^— ^^(g)ul 



39' (p(x)eV2u) + V {p{x)@^l'^vf^ 



(119) 



which can be written as 
dpu 



^ + div (pu (g) u + J5 (g) u) = pf - VP 

39* (p(x)eV2u) + V (p(a;)ei/2ti^) 



+ y9.p 



(120) 



4 Appendices 

4.1 Appendix 1 

The difference between the cases with and without an external potential i 

k' A;' [(A;'2 + k2)1/2 + ^/] ^ ^/ [(^/2 + ^2)l/2 _ ^/] 

m 2m(A;'2 + ^2)1/2 " 2m(A;'2 + ^2) V2 

^ (fc^2 + ^2)1/2 _ 

2m 

So we need to bound integrals of the form 
Bi := — exp <^ - 



IS 



2m 



Ck' )■ dk'. 



We have 



(fc'2 + ^2)1/2 _ ^/ 



< 



A:' + (A;'2 + ^2)1/2 

<2j^i if ^ 

K if < /e' < K. 
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Hence 



2mZie \Jo 



2m 



exp < - 



< 



< 



2mZj€ 



2 , 2 

K + K 



2m 

°° dk' 
-exp 



12 



2mZje 



- log 



k' y 2m 

1/2 



— I K 

m 



+ K 



^"^exp{l/2y2}dy^ 



Since = 2lmd^ and the term N^jiZie) is 0(1), we see that the error 
is 0(£log£). 

4.2 Appendix 2 

Here the difference between the cases with and without an external potential 
is 

k k-jk^ + K^/^ _ k+{e + KY^^ 

m 2m 2m ' 

so the integral we need to bound is 

f {k'^ + K^y/^ + k 



Bo :-- 



Zie 



k<0 



2m 



exp < -— Ck dk. 

2m 



With the change of variable k' = —k, we are led to the problem of finding 
an upper bound for 



^i^ Jk' 



+ ^2)1/2 _ y 



>o 2m exp^-^ + CV^d.'. 



f3k 



12 



But this reduces to the case of Appendix 1, since the quadratic (negative) 
term in the exponent eventually (and in fact very quickly, since m is so small) 
overcomes the linear (positive) one. 
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4.3 Appendix 3 



As in Appendix 1, the difference between the cases with or without an 
external potential is 



m 2m{kf + k2)1/2 



< 



2m{k'^ + k2)1/2 



2m 



so that we need to bound the integral 



k[>o ^rn 

,ikf + K^y/^-k[ 



A{k')d^k'. 



dk' 



2m 



k[>0 

k'k' 
2m 



which gives rise to the following four terms 



+ ) exp <^ - 



/3k' • k' 
2m 



-C-k' 



Zie 



+ 

+ 



*:i>0 



2m 



2m 



exp 



/2 



2m 



C'fe'i ^ dk[ 



(Zie)(2mZ2e) 

Ar 

(Zie)(2mZ3e) 



5i[M2(C') + M2(-C')] 
i?i[M2(C') + M2(-C')] 



Zie 



Using that, to zeroth order in C, M2 and Z^e are proportional, respec- 
tively, to m'^/2 and m~^^'^, we conclude that the last three terms are all of 
the same order as Bi, that is, 0{£log£). 

For the first of these terms, using the same estimate as in Appendix 1, 
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we find 



< 



< 



Am?Zie 



2m 



K j k-^ exp 





/2 



Pk^ 

2m 



exp < 

2m 2m 



Ck[ \ dk[ 



Ck[ \ dk[ 



+ K 



k[e^^{-^-Ck[\dk[ 



Am^Zie 



2m 



Zie 



W 



making it of a smaller order than the last three. Therefore, itself is of 
order 0{ilog£). 



4.4 Appendix 4 

The integral in (|69|) is bounded by 



2mZ€^ 



^(k')(i^k' 



lk[>0 



2mZ€^ .lk[>o V 2m 



^(-'))exp^-^-C-k'^dV 



2m 



2m 



2mZe^ 
2m 



Now, to zeroth order in C,, this reduces to 

4\/2/33/2mi/2(Zie) ^ 272/3^7w72(^ ^ 2V2^3/2mV2(Z3e) ^ ^ie2m 

Recalling that k? = 2imd^ and that the term Nx/{Zie) is 0(1), we find 
that the expression above is 0{l'^). 



28 



4.5 Appendix 5 

As in Appendix 2, the difference between the cases with or without an 
external potential is 



2 



m 2m 2m 

so that the integral to be bound in this is 



+ <^{x) exp <^ C • k M^'k 



2m J I 2m 

So performing the change of variable k' = — k, the integral becomes 



B. 



5 





1 




Jk[>0 


< 


k' k' 


2m 



2m 



+ ^(x) ) exp <! -^Kl^ + C • k' [> d^k' 
J [ 2m J 

which reduces to the case dealt in Appendix 3, since the quadratic term in 
the exponential quickly dominates the linear one. 

4.6 Appendix 6 

Replacing sums by integrals and using that Px{k) = Nxpxik) we find that 
the last two terms in (|102| ) can be expressed as 

^ ^/ \ ^x+eej f J Px+iejk'^ , I 

di^(x) — / exp < Ci,i>^ k} dk 

' ^ 'Zj{x + eej)eJk<o \ 2m '-+^^^ J 

f""';^ / expl- ^"-'^^^' -C^fe.feUfc 
Zj{x -ej)e Jk>^ [2m ^ '^^^ J 

We first notice that we can move the lower limit of integration in the 
second term above from k to zero, since the error involved in doing so is of 
order k, that is, 0(^3/2)_ -y^^ 

are then left with 

9jHx) y!^^l'7 , Mo{-Ci+,,.) + d,^x-£ej) JJ''-''; Mo(C^_,,.). 
Zj{x + £ej)e ^^^^ Zj{x — £j)e 
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Expanding Mq to zeroth order in 
'''^^''>Zjix + eej)e\2(3x+eeJ ^ ' '''' Zj{x - lj)e [20,^,^,) 



In the limit £ ^ this reduces to 

4.7 Appendix 7 

Here the difference is 

k' + (fc^^ + kY^^ k' _ jk'^ + K^/^ - k' 

2m m 2m 

and the integral to bound is 

Nx f + ^2)1/2 _ ^/ r f^y2 
Bt.= ^ ^ A;'exp<^-^^ Ck' \ dk' . 



ZjC Jf,i 2m y 2m 

Using the same estimate as in Appendix 1, we find 



which we conclude is 0(^^), since, to zeroth order in C, Mq is proportional 
to m^/'^. 

4.8 Appendix 8 

The integral to be bound here is 

:= — / A;exp<^-- (k'fdk. 

With the change of variables k' = —k this becomes 



AT, r k' + {k'^ + /.^)V^ ,^ r 
= -— / fcexp <^ Ck }dk , 



which has the same bound as that of Appendix 7. 
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